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Abstract The very weak solution of the Poisson equation with boundary data is 
defined by the method of transposition. The finite element solution with regularized 
boundary data converges with order 1/2 in convex domains but has a reduced conver¬ 
gence order in non-convex domains. As a remedy, a dual variant of the singular com¬ 
plement method is proposed. The error order of the convex case is retained. Numerical 
experiments confirm the theoretical results. 
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1 Introduction 

In this paper we consider the boundary value problem 

— Ay = / in P, y = u onT = dQ, (1-1) 

with right hand side / € and boundary data u € L^(r). We assume P C 

to be a bounded polygonal domain with boundary T. Such problems arise in optimal 
control when the Dirichlet boundary control is considered in T^(r) only, see for example 
the papers by Deckelnick, Gunther, and Hinze, [7], French and King, [8], and May, 
Rannacher, and Vexler, [11]. 
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For boundary data u £ F^(r) we cannot expect a weak solution y £ Therefore 

we define a very weak solution by the method of transposition which goes back at least 
to Lions and Magenes [10]: Find 

y £ L^(0) : (?/, Au)o = {u, dnv)T - (/, v)^ 'iv £V (1.2) 

with {w,v)g ■= JgWV denoting the L‘^{G) scalar product or an appropriate duality 
product. In our previous paper [1] we showed that the appropriate space V for the test 
functions is 

V := n with ffi(Q) -.= {v e Av g L^(n)}. ( 1 . 3 ) 

In particular it ensures dnV £ L^(F) for v G V such that the formulation (1.2) is well 
defined. We proved the existence of a unique solution y G L‘^{Q) for u G L^(F) and 
/ £ and that the solution is even in iL^/^(II). The method of transposition is 

used in different variants also in [8, 2, 4, 3, 7, 11]. 

Consider now the discretization of the boundary value problem. Let T/i be a family of 
quasi-uniform, conforming finite element meshes, and introduce the finite element spaces 

Yh = {vhGH\Q):vh\T^riyTGTh}, Yoh = Yh H Y^ = Yh\an. 

Since the boundary datum u is in general not contained in Yj^ we have to approximate 
it by L^(r)-projection or by quasi-interpolation. We showed in [1] that we can construct 
in this way a function with 

Ik - '“^lk-i/2(r) < 

As a side effect, the boundary datum is regularized since G Hence we can 

consider a regularized (weak) solution in Y^ := {u £ : u|r = u^}, 

y^GY!:-. iVy\Vvh = {f,v)n G (1.4) 

The finite element solution yh is now searched in := Y^ n Y^ and is defined in the 
classical way: hnd 

yh^y*h- {y'yh,y'vh)n = {f,Vh)n 'ivhGYoh. (1.5) 

The same discretization was derived previously by Berggren [2] from a different point of 
view. In [1] we showed that the discretization error estimate 

Ik “ yh\\L^(n) < Ch^ + lkllL2(r)^ 

holds for s = 1/2 if the domain is convex; this is a slight improvement of the result of 
Berggren. 

Let us now consider non-convex domains. Although the very weak solution y is also 
in 77^/^(n) the convergence order is reduced; the finite element method does not lead 
to the best approximation in L^(fl). In order to describe the result we assume for 
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simplicity that fl has only one corner with interior angle u € (tt, 27r). We proved in [1] 
the convergence order s G (0, A —^), where A := and showed by numerical experiments 
that the order of almost A — ^ is sharp. 

In this paper, we modify the discrete solution yh from (1.5) in order to retain the 
convergence order s = ^. In particular, we suggest to compute a function 

Zh Span{r“'’' sin(A0)}, 

where r, 6 are polar coordinates at the concave corner, such that the error estimate 

\\y - Zh\\L^(n) < + ll'«llL2(r)) 

can be shown. This method is a dual variant of the singular complement method intro¬ 
duced by Ciarlet and He [5]. Numerical experiments confirm the theoretical results. 

2 Analytical background and regularization 

As in the introduction, let be a domain with exactly one concave corner, and denote 
this interior angle by w G (vr, 27r). This corner is located at the origin of the coordinate 
system, and one boundary edge is contained in the positive xi-axis. It is well known 
that the weak solution of the boundary value problem 

— Av = g in Q, u = 0 on T = dQ, (2-1) 

with g G L^(n) is not contained in but in 

n n 0 Span{^(r) sin(A0)}, 

^ being a cut-off function, see for example the monograph of Grisvard [9]. This means 
that 

R-.= {Av.ve H^{Q) n 

is a closed subspace of L?‘{VL). It is shown in [9, Sect. 2.3] that 

^2(12) = i?0Span{p,}, (2.2) 

with the dual singular function 

Ps = r~^ sin(A6') + ps (2.3) 

where ps G is chosen such that the decomposition (2.2) is orthogonal for the 

L^(I2) inner product. Therefore, the dual singular function p^ is a solution of 

wGL'^in): {Av,w)=0 'iv £ H^{n) n (2.4) 

which proves the non-uniqueness of the solution of (2.4). This is the dual property to 
the non-existence of a solution of (2.1) in H^{Q) n see [9, Introduction]. 
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Due to (2.2) we can split any L^(D)-function into L^(D)-orthogonal parts. To this end 
denote by Ilij and the orthogonal projections on R and on Spanjps}, respectively, 
i.e., for g G L^(D), it is 5 = where 


^Ps9 = a{g)ps with a{g) 

nR5 = a- '^ps9- 


{9,Ps)n 

l|Ps|||2(Q) 


Since ps G L?‘{Q) there exists 


(t>s e n Hl{n ): -A0, = 


(2.5) 


see also Section 4 for more details on . For the moment we assume that ps and are 
explicitly known; hence the decomposition g = ling + a{g)ps can be computed once g 
is given. Computable approximations of ps and 4>s are discussed in Section 4. 

Now we come back to problem (1.2) and decompose its solution y in the form 


y = nRy + a{y)ps. ( 2 . 6 ) 

From the decomposition (2.2) we see that problem (1.2) is equivalent to 

{y,Ps)n = -iu,dn(j)s)r + {fAs)n, 

{y, Av)n = {u, dnv)r - (/, v)n Vu G F^(D) n Fd(D) 

and with the orthogonal splitting (2.6) to 

a{y) {Ps,Ps)n = -{u,dn(t)s)r + (/, </>s)o, 

{Uny, Av)n = {u,dnv)r - {f,v)n Vu G H^{Q) n 

The first equation directly yields a{y), namely 

^ ^ -{U,dn(j)s)r + {f,(psh 

{Ps,Ps)Q 

hence the projection of y on ps is known. It remains to find an approximation of Xl/jy. 

At this point we recall the regularization approach from [1] which we summarized 
already in the introduction. Let G be a regularized boundary datum such 

that we can define the regularized (weak) solution in := {v G Lf^(D) : u|r = 

y^^Y!:-. iVy\Vvh = {f,v)n e ( 2 . 8 ) 

In [1] we showed that the regularization error can be estimated by 

\\y - /||l2(o) < c\\u - u^\\H-'>{r) 

where s = ^ if D is convex and s G [0, A — ^) if D is non-convex, that means the 
regularization error is in general bigger in the non-convex case. With the next lemma 
we show that n/j(y — y^) is not affected by non-convex corners. 
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Lemma 1. If the domain It is non-convex, the estimate 

l|n/j(y - /)||L 2 (f^) < C\\u - 'u'^||^-i/ 2 (r) 

holds. 

Proof. Recall V = n Hq{Q) from (1.3). From (2.8) and the Green formula, we 

have for any v € V 

{f,v)n = (Vy'", Vn)o = -{y’^,Av)n + {y'^,dnv)r. 

Note that n G R is sufficient, see [6, Lemma 3.4]. Subtracting this expression from the 
very weak formulation (1.2), we get 

{y-y'^,Av)n = {u-u’^,dnv)r Vn G R. 

Restricting this identity to n G H‘^{Q) n Hq{Q), we have 

{Uniy-y^),Av)n = {u-u^,dnv)r Vn G iL^(O) n ifg (12). (2.9) 

Now for any z £ R, we let Vz G n Hq{PI) be the unique solution of 

Avz = z, (2.10) 

that satisfies 


\\dnVz\\H^/2^Y) < cWVzWh^Q) < cW^L^iU)- 
Since for any g G Lp‘{^) the equality 

(^R{y - y^),g)n = (nij(y - y'"),URg)n = (y - y^, 'nR9)i 

holds we get with (2.9)-(2.11) 


( 2 . 11 ) 


W^Riy-y'^)\\L^{n) = sup 


{y-y^,z)n {u-u^,dnVz)r 


‘U II II — sup .. .. 

,z^0 |R||l2(q) zeR,z^0 IRIlL2(f2) 


zeR. 


< Wu-u'’ 


ll/fi/2(r) 


|j^-i/ 2 (r) sup 

z£R,z^o W^UL^in) 

which is the estimate to be proved. 


< cllu — u 


h I 


l//-i/2(r) 


□ 


3 Discretization by standard finite elements 

Recall from the introduction the finite element spaces 

Yh = {vh£H\n):vh\TeViyT£Th}, Yoh = Yhn H^{n), Y^ = Yh\dn, 
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defined on a family Th of quasi-uniform, conforming finite element meshes. Assume that 
the regularized boundary datum is contained in such that the estimates 

||w^||L2(r) < c||n||j;^2(r), (3.1) 

11'^ ~llH-i/2 (r) ^Ch / , (3-2) 

hold. It is proved in [1] that this can be accomplished by using the L^(r)-projection or 
by quasi-interpolation. A consequence of Lemma 1 is the estimate 

l|nij(2/ - /)|U2(o) < Ch^/^u\\L2^r) (3.3) 

in the case of a non-convex domain Q. (In the case of a convex domain the operator 
is the identity, and the corresponding error estimates were already proven in [1].) 

As already done in the introduction, define further the finite element solution yh £ 
:= y*'^ n Yh via 


Vh £ y*h ■ (Vyfe,Vu/,)o = {f,Vh)n Vu/* £ Y^h- (3.4) 

We proved in [1] that 

\\y — yh\\L^{n) < + ll^^llL 2 (r)^ (3.5) 

holds for s = ^ if the domain is convex but only s £ (0, A — ^) in the non-convex case. 
In the next lemma we show that nj^(y — yh) is not affected by the non-convex corners. 

Lemma 2. For non-convex domains II the discretization error estimate 
l|nij(y — yh)\\L^{n) < + ll^llL2(r)) 

holds. 

Proof. By the triangle inequality we have 

||n/j(y - yh)\\L^{n) < l|nij(y - y’^)\\L^(n) + ||n/j(y'' - yh)\\L^{n)- (3-6) 

The first term is estimated in (3.3). For the second term we first notice that £ 

Hq{Q,) satisfies the Galerkin orthogonality 

{V{y^ -yh),Vvh)n = 0 £ yoh, (3.7) 

see (1.4) and (1.5). With that, we estimate ||n/j(y^ — yh)\\L'^(n) by a similar arguments 
as ||nj:j(y — y^)||L 2 (Q) in the proof of Lemma 1. Recall from (2.10) and (2.11) that 
Vz £ iL^(fl) n Hq{Q) is the weak solution of Avz = z G R. It can be approximated by 
the Lagrange interpolant Rvz satisfying 

||V(U^ - IhVz)\\L^{Q.) < C.h\\Vz\\H^{Q.) < ch\\A\L^{Q.)- 
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We get 


W^niy’" - yh)\\mQ) = sup 


{y'^ - yh,z)i 


= sup 


(V(y'* - yh)yv^)n 


^ '^Jr 11 11 ^ 

z€R,z^O ll^llL2(r2) z&R,zj^O 


"llL2(n) 


= sup 

z£R,zj^O 


(V(y^ - yh),^{vz - hvz))i 


'llL2(n) 


< ch\\V{y^ - yh)\\L2{n)- 


(3.8) 


In order to bound ||V(y^ — yh)\\L^{n) by the data we consider a lifting G 

; 8 

= <;: r T ^( 3 - 9 ) 


defined by the nodal values as follows: 

u^{x), for all nodes x G F, 


0 


for all nodes x G n. 


The homogenized solution yg = G Hq{Q) satisfies 

(v7,Vu)o = if,v)n - {V{BhU^),Vv)n Vu G 
By taking v = y^ we see that 

l7yolli2(f^) < Il/I7-I(o)l|yol7i(f^) + \\"^i^hu'^)\\L^{n)\\"^yo\\L^{n)- 

Using the Poincare inequality we obtain 

INyohnn) < c||/l7-i(c) + ||v(b,u")|U 2 (^), ( 3 . 10 ) 


and with the Cea lemma 

l|V(7 - yh)\\LHn) < l|V(7 - BhU^)U2^n) = I|v7||l2(o) 

The remaining term ||V(i?/itt^)||^ 2 (Q) is estimated by using the inverse inequality 


\\'V{Bhu’^)\\L^(T) < ch 

for £■ C r n T, T G Th, which can be proved by standard scaling arguments, to get 

||V(B;,u'‘)||i2(o) < ch-^/V\\L^(r)- (3.11) 


Hence we proved 

17 ( 2 /^ — yh)\\L^(n) < c\\f\\H-i{n) + ch ^/^||u^||i2(r). 

With (3.6), (3.3), (3.8), the previous inequality, and (3.1) we finish the proof. □ 

With (2.6) we can immediately conclude the following result. 
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Corollary 3. Let Q be a non-convex domain and let yh S be the solution of (3.4), 
then the discretization error estimate 

\\y - (nijyh + a{y)ps)\\L 2 (n) < Ch^^’^ + \\u\\L^{r)) 

holds, reminding that ps and a{y) are given by (2.3) and (2.7), respectively. 

Hence the positive result is that TlRyh + a{y)ps is a better approximation of y than y^. 
The problem is that Ps and (pg are used explicitly, and in practice they are not known. 
A remedy of this drawback is the aim of the next section. 


4 Approximate singular functions 

Following [5], we approximate ps from (2.3) by 

V^s=Ph-f'h + r~^ sin(A6»), r^ = (r~^ sm{Xefj , 

p*hGYoh- {Ypl,Vvh)n = {'^rh,Vvh)n V?;/i G Tq/i, 
with Bfi from (3.9). The function (pg from (2.5) admits the splitting 

(pg = ^ + fdr^sm{X6), 


(4.1) 


(4.2) 


with (p G and /3 = tt ^||Ps|li 2 (f^p see again [5]. It is approximated by 

(pg = (Pl-/3hSh +Phr^sm{X9), st = Bh (r^ sm{Xe)'^ , l^h = ^WPgWh^^p 

(pl^Yoh- iY(pl,Vvh)n = {Ps,Vh)n + Ph{Y'sh,Yvh)n Vu/i G Tb/u 

that means, (p is approximated by cph = (p\ — PhSh G Y^. The approximation errors are 
bounded by 

WPs - Pgh^in) < < ch, (4.4) 

1/3 - Ph\ < < ch, (4.5) 

\\(ps - fffWi,Cl < ch, (4.6) 

see [5, Lemmas 3.1-3.3], where (4.5) and (4.6) imply 

\\f - ^h\\i,n < ch. (4.7) 

At the end of Section 3 we saw that + a{y)ps is a better approximation of y 

than y/j. Since this function is not computable we approximate it by 


Zh = y/i + ahPg, 


(4.8) 


with 


^Ryh = yh-ihPs^ ih 
and a suitable approximation ah of 


{yh,Psh 

mUn) 


(4.9) 


a{y) 


-{u,dn(j)s)r + (/, 4 's)n 
{Ps,Ps)n 


from (2.7). To this end we write the problematic term by using (4.2) as 

{u, dn(j)s)r = {u, dn^)r + I3{u, (9„(r^ sin(A 6 i)))r. 

and replace the term {u,dn^)r by {u^,dn^)r- Since ^ belongs to and is the 

trace of BhU^, we get by using the Green formula 


{u^,dn^)r = {Bhu\A^)n + {VBhu\V^)n 

= -{BhU^,Ps)n + (VBhU^, V(^)o (4.10) 


as A(p = A4>s = —ps- With all these notations and results, we define 

{BhU^,Ps)n - {VBhU^,V^h)n - /3h{u,dn{r^ sin(A 6 »)))r + {fA^s)n 


Oih = 


[ph^ph\2 


(4.11) 


sJn 


Note that ah can be computed explicitly and therefore Zh as well. 
Let us estimate the approximation errors made. 


Lemma 4. Let Ll be a non-convex domain and let Uh S Y^h be the solution of (3.4). 
Then the error estimates 

WB^Rdh — n^yh||L 2 (o) < ch (||/||h-i(o) + ||tt||L 2 (r)) ) (4-12) 

|a(y) - oih\ < + ||w||i, 2 (r)) (4.13) 

hold. 

Proof. With the definitions of Ilij and 11^, with 7 := {yh,Ps)n/\\Ps\\\ 2 (^Qy and by using 
the triangle inequality we have 

\\URyh - n^2//i||L2(o) = WlPs - IhPsh^in) 

< I7 - 7 h| WPsh^n) + I7I WPs -Psh^n) 


We write 


7 - 7h 


{yh,Psh _ {yh,Psh 

117*^1112(0) 117*51112(0) 


{yh,Ps -Psh 

WP^Wt^in) 

{yh,Ps -Psh 

117*51112(0) 


+ {yh,Ps)n 


1 



\\P 


h\\2 

5 IIl2(o) 


+ {yh,Psh 


{Ps +Ps,Ps -Ps)q 

ll7'5lli2(o)ll7'5 lli2(o) 
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and by the Cauchy-Schwarz inequality and (4.4) we get 

|7 -7h| < chWvhWL^i^n)- 

We have used that and can be treated as constants due to the 

dehnition of Ps and due to (4.4). We conclude with I 7 I < cWyhWi'^^n), and (4.4) that 

W^RUh - n^y/i||i,2(o) < chWuhWL^^n)- (4-14) 

In view of the finite element error estimate (3.5) and the standard a priori estimate for 
the very weak solution, 


\\y\\L^{n) < c (||/||j^-i(Q) -I- ||tt||i2(-r)), 
see Lemma 2.3 of [1], we have 

l|y/i||L2(r2) < l|y||L2(n) + \\y - yh\\L‘^(n) < C (||/||_H-i(f2) + ll“llL2(r)) • 

This estimate together with (4.14) proves (4.12). 

The proof of the estimate (4.13) is based on writing the problematic term in the 
definition of a{y) without approximation as 

{u, dn(l)s)r = {u, dn4>)r + /3(tt, dn{r^ sin(A6»)))r 

= {u- u’^,dn^)r + {u^, dn^)r + /3(u, dn{r^ sin(A6»)))r 
= {u- u’^,dn4>)r - {Bhu’",ps)n + {VBhU^, V^)n + (3{u, dn{r^ sin(A6»)))r 

where we used (4.10) in the last step. Consequently, we showed that 

a(y) -ah= n—li^-( - {u - u^,dn^)r + {BhU^,ps - p^s)^ - {VBhU^,V{^ - ^h))n 

\\Ps\\l 2 (^Q) ^ 

-{13- Ph) {u, dn{r^ sm{Xe)))r + (/, 4>s - (ps)nj ■ 


To prove (4.13), in view of (4.4), (4.5), and (4.6) it remains to show that 

{u-u^,dn^)r < c/i^/^||u||2,2(r), 
{Bhv!^,Ps-p^s)n < c/i^/2||u||i2(r), 

{VBhu\V{^-^^))n <chV2||^||^2(r). 


The first estimate follows from the estimate (3.2) and the fact that (j) belongs to 
The second one follows from the Cauchy-Schwarz inequality and the estimates (3.11) 
and (4.4). Similarly, the third estimate follows from the Cauchy-Schwarz inequality and 
the estimates (3.11) and (4.7). □ 
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Corollary 5. Let LI be a non-convex domain and let yh G he the solution of (3.4) 
and let Zh be derived by (4.8), (4.9), and (4.11), then the discretization error estimate 

\\y — ^/iIIl 2 (q) < + ll'“llL 2 (r)) 

holds. 

Proof. The main ingredients of the proof were already derived. Indeed, it is 

\\y- Zhh^n) = W^ny + a{y)ps - IlRyh - ahPs h^n) 

< linijy - - n^yhllL2(o)+ 

|a(y) - ah\ WPsWmn) + l«/i| WPs-Psh^n)- 

The hrst three terms can be estimated by using Lemmas 2 and 4. So it remains to 
treat the fourth term. To bound \ah\ we use the triangle inequality 

\o!h\ < \ah - a{y)\ + \a{y)\. 

For the first term we use (4.13), while for the second term we use (2.7) reminding that 
(ps belongs to with some e > 0. Altogether we have 

\oih\ < C (||/||h-i(o) + ||^f||L2(r)) 

and conclude by using (4.4). □ 

Before we describe the numerical experiments, let us summarize the algorithm. 

1 . Compute the finite element solution 

yh G : (Vy/i, Vuh)o = (/, Ufe)^ G Tbh 

where = {vh G : u/i|r = u^}-, compare (1.5), with G Y)f being an 
approximation of the boundary datum u satisfying (3.1) and (3.2). 

2 . Compute the approximate singular functions: 
rh = Bh (r-^sm{Xe)^ , 

Ph G Yoh ■■ (Vpl, Vvh)n = C^rh,Vvh)n 'dvh G Tq/i, 

Ph=Ph- i"h, 

Ph = ^||pft + r-^sin(A0)|||2(t^), 

Sh = Bh [r^ sin(A6l)^ , 

(pl^Yoh- {V(p*h,Vvh)n = {ph + r~^sm{\e),vh)n +PhiY'sh,Y'vh)n Vu/* G Tbh, 
Ph — Ph ~ PhSh) 
compare (4.1) and (4.3). 
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3. Compute 


Ih = 


Oih = 


{yh,Ps)n 


with Pg = pfi + r ^sin(A0), 


{Ps,Psh 

{BhU^,Ps)n - C^BhU^,V^h)n - Phiu,dn{r^ sin(A6i)))r + (/, (j)s)n 

{Ps,Ps)l 


Q;/i T/i) 

Zh = yh + ShPh, 

compare (4.9) and (4.11). According to (4.8), the numerical solution is 

Zh = Zh + 6hr~^sm{Xe). 

Note that all integrals with and r~^ must be computed with care. 


5 Numerical experiments 

This section is devoted to the numerical verification of our theoretical results. For 
that purpose we present examples with known solution. Furthermore, to examine the 
influence of the corner singularities, we consider several polygonal domain depending 
on an interior angle w G (0, 27r). The computational domains are defined by 

Xlu) ■= ( — 1; 1)^ n {x G : (r(x),0(x)} G (0, V2] x [0, w]}, (5.1) 

where r and 9 stand for the polar coordinates located at the origin. The boundary of 

is denoted by F^. We solve the problem 

—Ay = 0 in y = u on F, (5-2) 

numerically by using the proposed dual singular function method. The boundary datum 

u is chosen as follows 

u := sin(-O.49990) on F^. 

This function belongs to LP{T) for every p < 2.0004. The exact solution of our problem 
is simply 

y = sin(-0.49996»), 

since y is harmonic. 

The quasi-uniform finite element meshes for the calculations are generated by using a 
newest vertex bisection algorithm. The discretization errors for different mesh sizes and 
the corresponding experimental orders of convergence are given in Table 1 for different 
interior angles uj = 270° and uj = 355°. We see that the numerical results confirm the 
expected convergence rate 1/2. 

We emphasize that the quadrature formula for the numerical integration of the integral 

{u,dn{r^ sin(A6')))r 
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mesh size h 


eoc 

0.25000 

0.58725 


0.12500 

0.42338 

0.47201 

0.06250 

0.30318 

0.48177 

0.03125 

0.21606 

0.48870 

0.01562 

0.15352 

0.49302 

0.00781 

0.10888 

0.49572 

0.00390 

0.07712 

0.49742 


mesh size h 


eoc 

0.25000 

1.02069 


0.12500 

0.83402 

0.29139 

0.06250 

0.58964 

0.50025 

0.03125 

0.41696 

0.49991 

0.01562 

0.29506 

0.49890 

0.00781 

0.20903 

0.49725 

0.00390 

0.14836 

0.49462 


Table 1: Discretization errors eh = y — Zh ioi u = 37r/2 (left) and cj = 3557r/180 (right) 


has to be adapted in order to get a sufficiently good approximation. Otherwise, the error 
due to the quadrature formula dominates the overall error. In our implementation, we 
chose for the numerical integration a graded mesh on the boundary {He r\j ^ if the 
distance te oi the boundary edge E satisfies 0 < < i? with R being the radius of 

the rehnement zone and y being the rehnement parameter, and for te = 0) 

combined with a one-point Gauss quadrature rule on each element. Furthermore, the 
grading parameter y is chosen such that 

jjL < 2 'k/u} — 1 , 

which seems to be the correct grading to achieve a convergence order of 1/2. For the 
results presented in Table 1 we used i? = 0.1 and y = 2 t:/ uj — 1. 
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